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I. INTRODUCTION 


As the uniquely distinguished feature of quantum mechanics [T], and, as the primary resource of quantum informa¬ 
tion processing (QIP) [5] a deeper understanding of entanglement carry both foundational [5] and practical values, 
such as the application to quantum teleportation [4]. The first major task of finding mathematically sound [5] def¬ 
initions to quantify the existence of bipartite and multipartite entanglement for discrete variable quantum systems 
such as ‘concurrence’ [5] and continous variable quantum systems such as ‘negativity’ [7] was undertaken in the 90s 
From physical considerations a system where quantum information processing is carried out is always under 
the influence of its environments, whose effects can be represented as noises of different nature which could be detri¬ 
mental to QIP. Environment-induced quantum decoherence and disentanglement are two vital processes which need 
be accounted for and understood well to enable one to find ways to mitigate or manage them. This second major 
task is most effectively undertaken with the methods of quantum open systems m based on various aspects of the 
quantum Brownian models. Environments when suitably prepared or attuned to can also assist in maintaining or 
even generating entanglement mun]. Serious studies of entanglement between two qubits [TMTS] and continous 
variables UMl] appeared in the last decade, Our studies on how entanglement in a quantum system is affected by 
its environment follow this call, relying mainly on functional methods [23I126j . and in the case of Gaussian systems 
[?7l - E5] . the covariance matrix theory . 

In this paper we consider the entanglement between two harmonic oscillators, call detectors [? ], which make up 
our system, interacting with a quantum scalar field acting as its environment [ 301131 ] with special attention paid to 
two factors: distance and coupling dependence of the entanglement dynamics. 

One important factor which has been largely ignored or grossed over is how entanglement between two qubits 
(representing discrete variables) or harmonic oscillators (representing continuous variables) depends on their spatial 
separation. This is one of the two primary aims of our prior studies mill]. Distance dependence is recognized as 
an essential factor in relativistic quantum information |32j such as entering in quantum teleportaion considerations 
[35] . Even with stationary detectors this is a relevant factor, because the two detectors are linked by a field, and 
information propagates through the field as their medium, whereby entanglement can be induced, even generated. 
Many prior work (IMSH] assumes zero spatial separation between the two detectors and report on the (common) 
environment-induced entanglement dynamics - such as generation, death and revival. One may casually think that 
the distance effect is not so important because at most it invokes a rather weak retardation effect from the field 
environment, or a small contribution from the bath correlations. But a moment’s thought calls into question the 
following: If there is no direct coupling, as assumed, how could two detectors get entangled via a field which acts only 
at the same point where the two detectors are located? Both detectors would read off the same local value of the 
field, with no dynamics of the field involved. The environment would affect each detector’s entanglement dynamics 
but would not enter into influencing the combined two-detector system over and beyond the simple sum of each. 
Equivalently this amounts to a contact potential with instantaneous action of the field. Since this is the only means 
of interaction between the two detectors, the zero separation is clipping the wings of the field and thus the setup has 
some intrinsic deficiencies. 
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For an oscillator bath, one may argue that even without spatial correlation there should still be time-time auto¬ 
correlation as influenced by the environment (via the self-induced non-Markovian dissipation kernel). But the fact 
is, each detector would get the same influence from the environment, as it interacts with the environment in the 
same way as the other does. Without spatial separation or without direct coupling this would be like two separate 
detectors each conducting its own business albeit interacting with the same environment. The environment certainly 
influences each detector, but it does not contribute to any additional inter-oscillator entanglement beyond which it 
imparts to each detector. For these reasons, conclusions drawn from calculations making this assumption can only 
be valid for a very short time scale, reflecting a transient behavior which cannot be extended to late times. Thus the 
meaning conveyed by the word ‘common’ in “entanglement induced by a common bath” implying that there is added 
inter-detector entanglement due to the fact that the two detectors share a common bath could be misleading. The 
proper treatment even for the zero separation case would have to begin with a finite separation, give full dynamics 
to the field, allow it to partake of the dynamics of the combined two-detector system, then take the small separation 
limit. Indeed in our present calculation we can pinpoint and expose this shadowy feature and show the cause of it 
explicitly. This is one exemplary observation which calls for a careful reconsideration of prior results based on the 
zero-separation condition. 

The other important factor is allowing for direct coupling between the two detectors. This enlarges the scope from 
cases with no direct coupling, such as in neutral atoms, to include a wider range of inter-detector interactions. The 
type of coupling we treat here in detail is of the quadratic type between the detector’s internal degrees of freedom 
(denoted as y here, or Q in our earlier papers pTl [5^ 155] 1. which is of the same form where the HPZ-type of master 
equation [5SJ |3S] was earlier derived for two coupled oscillators m- By keeping the whole system Gaussian, we are 
able to produce exact analytical solutions valid even for strong couplings. Our present work is a direct descendent 
from and a generalization of our earlier work m on the temporal and spatial dependence of entanglement with the 
same setup, but now with non-zero inter-detector coupling. This allows us to see both a) the competition between 
the direct coupling in protecting entanglement and the environmental influence of corrupting entanglement and b) 
the direct coupling effect versus the indirect or induced effect of the field on the entanglement between the detectors, 
the latter came alive from a finite separation. 

We add that these considerations are not just adding details to established prior work, they surely do, but more 
significantly they enable us to identify some largely unnoticed blind spots in earlier treatments, to rectify possible 
misleading conclusions affected by them, and to provide a more reliable foundation for further developments. Including 
these two essential factors - distance and coupling - in the consideration of the entanglement dynamics of this 
workhorse model for continuous variable quantum systems and its extension could provide a more complete, more 
correct and more accurate analysis for many problems it can tackel with, such as expounding some quantum foundation 
issues (e.g., [38]), quantum information and teleportation applications (e.g., jSS] l39|), understanding macroscopic 
quantum phenomena (e.g., HOI HD) and exploring the interesting new field of quantum thermodynamics (e.g., lu¬ 
lls]). 


Below we first provide the background by describing some representative prior work in their respective ranges and 
the claimed results. We then provide a qualitative description of the main features of this system, present a brief 
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summary of our findings and note the main differences from prior work. 


A. Prior Works 

To provide the backdrop for this investigation we tabulate what has been done before in hve oft-cited papers: 
Columns in Table |l] indicate the authors corresponding to References [T8H^ respectively; Rows A-H indicate whether 
distance and coupling between the detectors are considered, properties of the bath, etc. The footnotes beneath the 
Table signify some special conditions related to a specihc feature. A quick glance of Items A and B indicate that the 
two papers closest to our concern are that of [20] and [21], factoring in much overlap of contents in [22] with [2T]. 

B. Direct and induced interactions, Symmetry Considerations, Separation and Coupling Dependences 

We assume that two detectors stay at rest at a hnite distance t apart in a common bath modeled by a massless 
scalar field. The internal degree of freedom y of each detector is described by a harmonic oscillator, which is directly 
coupled to the internal degrees of freedom of the other oscillator and interacts with the environmental scalar field, both 
bilinearly. To assess the validity of the results reported in prior and the present work it is necessary to understand 
how different interactions in the system affect the dynamics and how its behavior varies in the different parameter 
regimes of interest, as we now describe. 

I.l. Induced interaction and non-Markovian Dynamics The two types of interactions between the internal degrees 
of freedom of the two detectors, a direct inter-oscillator coupling, and an indirect, nonlocal, retardation interaction, 
mediated by the common environment field determines the dynamics of the system quite differently. We can see 
why the indirect interaction produces non-Markovian dynamics as such: the interaction of one oscillator with the 
field produces some disturbance in the field (or signal) which takes a finite amount of time (follows the lightcone) to 
propagate to the other oscillator located at some distance and modify its motion. The second oscillator does the same 
in emitting its signal while receiving the signal from the first oscillator, but at a later time. The history-dependent 
nature is what gives the induced interaction a non-Markovian character. This is explained and well-illustrated in the 
results of our prior work |21] . 

We may add a few words on the difference in the role of the environment made up of a scalar field (e.g., [SO] |3T] 
) versus fV-harmonic oscillators (e.g., [531 |M]. To see this relation explicitly, as e.g., in [^, where a scalar field is 
explicitly expressed in terms of parametric oscillator modes. But there are subtle differences: We mention just two 
here. First, the spectral density of a field is fixed, and takes on a different form for different dimensions. However, 
one can specify the form of the spectral density function for an A-oscillator bath. Second, how to deal with a 
divergent contribution from the high frequency end of the spectrum. There are well-justified ways to obtain a well- 
dehned integral over the whole frequency domain in the case of a quantum field bath, by way of regularization or 
renormalization methods. For an 7V-oscillator bath the spectral density is often crafted to take on a form that its value 
tapers off to zero when the frequency is greater than a certain high frequency cutoff scale. As a consequence, in the 
time domain the dissipation kernel will subside at a rate of the inverse cutoff scale. For a very high frequency cutoff 
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Comparison Table*-^^ 

Liu-Goan |18| 

Horhammer- 

Biittner 

m 

Paz-Roncaglia 

Lin-Hu 1 ^ 

Zell-Queisser- 

Klesse 

m 

A. Spatial separation 

N 

N 

N 

Y 

Y 

B. Direct coupling 
between detectors 

N 

N 

yC) 

N 

N 

C. Bath type 

A-oscillator 

A-oscillator^®^ 

A-oscillator 

scalar held 

Extended 

A-oscillator*'"'^ 

D Time regime 

Transient 

Transient 

Transient & late 
time 

Transient & late 
time 

Transient & late 
time 

E. Equation of 
motion 

Master Eq. with 
RWA-I-Markov 
approximation 

Exact master 
equation 

Exact master 
equation 

Quantum 
Langevin Eq. 

Quantum 
Langevin Eq. 

F. Initial state 

Two-mode 
squeezed vacuum 

Two-mode 
squeezed state 

Two-mode 
squeezed state 

Two-mode 
squeezed state 

ground state 

G. Computation 
approach 

Numerical 

Numerical 

Analytical & 
numerical 

Analytical 

Numerical 

H. Critical 
parameter 

None 

None 

Sc - ^ 

V 2moj 

d ~ 

wAi 

dc ~ A-i 


TABLE I: Comparison of five oft-cited papers. 

1. Symbols: m, w and 7 denote the mass, natural frequency, and damping constant of both oscillators, respectively. 
The parameter i describes the spatial separation, and cr represents the coupling strength between the two 
oscillators. 

2. A direct coupling constant is present in the formalism but was renormalized to zero in their calculation. 

3. The authors assume that the dynamics of the relative variable comes about as the solution of a Born-Markovian 
master equation for a damped oscillator while in actuality the relative variable really behaves like a free oscillator, 
a fact from symmetry considerations. 

4. Spatial dependence is inserted in the interaction term of the system oscillators and the bath oscillators, so in a 
sense the configuration is very similar to that in m 

5. The critical temperature is found for the coherent state with direct coupling strength renormalized to zero so 
it is independent of the squeezing parameters and the direct coupling. We will treat thermal entanglement in a 
sequel paper |45) . 

6 . The dimensionless parameter Ai is related to the resolution At of the detector by Ai = — In a; At — The 
critical separation is independent of the coupling constant only because it is obtained in the vanishingly weak 
coupling limit. 


this falloff is very steep, but for not so high a frequency cutoff, reflective of some ‘softer’ physical environment, the 
dissipation kernel would contain a feature which resembles a self-induced time-nonlocal interaction of each oscillator 
with itself. This feature is usually ignored under the assumption that fewer high-frequency modes in the environment 
take part in the interaction with the subsystems, from unspecified physical considerations. 

1.2. Symmetry and Dynamical Considerations One can regroup the two internal degrees of freedom of the two 
detectors into a symmetric (or center of mass) variable and an anti-symmetric (relative or difference) variable and 
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use them to describe the dynamics of the combined (two-oscillator) system. The dynamics of these two variables 
(interchangeably may be called ‘modes’) are very different and dictate very different entanglement behavior. (An 
explicit illustration of how symmetry determines the entanglement behavior is given in m for a system of two qubits 
interacting with two cavity modes.) In a A^-oscillator bath the relative mode is found to be freely evolving, totally 
decoupled from the bath. The decoupled relative mode dynamics is the root cause for the fact that the late-time 
entanglement of the combined system oscillates in time and shows dependence on the initial condition of the joint 
system, even though it is acted upon by a bath whose damping function one expects to wash out the sensitive 
dependence of the initial conditions. This entanglement behavior is unavoidable in the dynamics of two oscillators at 
zero-separation (attempts to remove this behavior by adding some damping by hand is unwarranted and unphysical). 
However, with finite separation between the two detectors, field-induced effects set in to provide a natural damping 
effect on the relative or anti-symmetric mode dynamics. From a Taylor expansion of the retardation terms we see 
immediately that the entanglement dynamics of the two oscillators at the same location shows up at very short times 
as transients. Thus our earlier comment that results obtained from the zero-separation setup cannot be taken to 
represent the long time behavior. Including higher-order terms in this expansion shows that the relative mode decays 
with time at a slow rate. For sufficiently long time after the joint system has reached full relaxation, the entanglement 
approaches a constant and becomes independent of the initial conditions. 

1.3. Oscillator-Bath Coupling and Instability For non-Markvian dynamics in the current configuration, strong 
oscillator-environment coupling will not result in an overdamped motion as it usually does for Markovian dynamics. 
Instead it tends to induce instability. This in a sense is similar to parametric oscillation: the disturbance from one 
oscillator can have a secular amplihcation effect on the motion of the other oscillator. 

More details can be found in the calculations described in the following sections. 

C. Temporal, spatial and coupling strength- marked regimes 

In the following, we assume the two oscillators have the same natural frequency oj and damping constant 7 . The 
parameter £ describes the spatial separation, and a represents the coupling strength between the two oscillators. 

• Short/long Time. The dividing time scale is the relaxation time 7 “^. We are mostly interested in the 

entanglement behavior of the two coupled oscillator system at late times, after it has fully relaxed or reached 
equilibrium with the environment. We assume that the oscillators have exchanged their (non-Markovian) mutual 
influence via the field many times, which implies that 7 “^ > i. For a small factor, the non-Markovian 

effects decays rapidly after a few exchanges. 

• Short/large Separation. When the distance between the two oscillators is shorter than 7 /w^, either the 
perturbation calculation is unreliable or the motion of the joint system becomes unstable. The instability is 
a feature of non-Markovianity different from the simplifying (but often lack of justification) assumption of a 
Markovian process leading to overdamped motion of the system. Another interesting scale is 7 /cr where a change 
in the entanglement measure occurs. It is the value where the direct interaction (inter-oscillator coupling) is 
comparable in strength with the indirect (field-induced) non-Markovian effect. When the separation between 
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the oscillators is larger than this scale, the direct coupling between the oscillators dominates the dynamics. On 
the other hand, when the separation is shorter than that scale, the non-Markovian (field-induced) effect is the 
determining factor in the motion of the joint system. Other than these two length scales, we refer to wf 1 as 
large separation and a;£ <C 1 as small separation. 

• Weak/strong Coupling. For the coupling between the bath and the oscillators, we only consider the case 
y/w < 1, that is, weak coupling. When 7 is comparable or greater than w, instability of the motion tends to 
occur unless 7 /w. This range is less interesting because it corresponds to the very large separation cases. 

For direct couplings between the oscillators, we only require that the resonance frequencies of the normal modes 
be real. This implies cr < — 7 ^, which is not a stringent constraint. Thus a can be comparable with up' and is 

considered strong. Thus we may refer to cr < 27 /^ as weak coupling because the direct coupling plays a minor 
role compared with the non-Markovian indirect interaction. In certain cases, we may encounter a < 70 ; < 27 /t!. 
It corresponds to vanishing inter-oscillator coupling. This is equivalent to the cases, in the table, where the 
inter-oscillator coupling is not present. 


D. Late time entanglement behavior 


For finite separation between the two detectors, due to the field-mediated (non-Markovian) mutual influence inter¬ 
action, the effective damping constant and the effective oscillating frequency are separation-dependent and do not take 
on the same values between the two (symmetric vs antisymmetric) normal modes. Being at late times, the motion 
of either normal mode is fully relaxed or has reached equilibrium. Since the strength of mutual influence depends on 
the separation between the two oscillators, and there is competition between the direct inter-oscillator coupling and 
the indirect, environment-mediated interaction, the late time entanglement between two oscillators can have quite 
distinct behaviors, depending on the ratio c = — ^ I, as follows: 

• When c > 1: this is the regime where the direct coupling wins over the indirect interaction between the 
oscillators. We find that the late time entanglement can improve with 


— stronger inter-oscillator coupling, 

— larger oscillator separation, and 

— weaker oscillator-bath interaction. 


• When c < I: the effects of direct coupling is subdominant, compared with the field-induced interaction. We 
find that the late-time entanglement is favored with 

— weaker inter-oscillator coupling, 

— shorter oscillator separation, and 

— stronger oscillator-bath interaction. 

Alternatively, heuristically speaking, for sufficiently strong inter-oscillator coupling, their late time entanglement can 
survive even when they are very far apart. The late time entanglement worsens if the two oscillators are originally 



placed closer. When their separation is less than a critical value £>, as shown in (124) or (130), the oscillators will 


instead become separable at late time. It means that direct coupling is not enough to overcome the effect of both 
the field-induced non-Markovian interaction effect and the local vacuum fluctuations of the field. When we place 


the oscillators at a separation even shorter than another critical value £< as in (138), the entanglement between the 


oscillators can again be sustained at late time. Now that the non-Markovian processes have dominated over those 
from direct coupling the results can be less intuitive and harder to interpret. If we set up oscillators at shorter and 
shorter separation, we find that their motion eventually become unstable. Quantitative results are contained in the 
last section. 

In short we see that direct coupling can play a key role in maintaining entanglement between the subsystems 
over finite distance. The type of inter-oscillator coupling we used may be idealized in the sense that its strength 
remains constant over large separation. Even when the inter-oscillator coupling is separation-dependent and say 
falls off relatively slowly with larger separation, the arguments above should still be valid over the region where the 
inter-oscillator coupling is much stronger than the mutual-influence interaction between the oscillators. 

The organization of this paper is as follows: In Sec. 2 we discuss the dynamics of two coupled quantum oscillators 
at finite separation as the system in a common quantum field environment, paying particular attention to their late¬ 
time evolutions. In Sec. 3, we construct the covariance matrix from the uncertainties or the cross-correlation of the 
canonical variables of these two oscillators. This enables us to compute the entanglement measure. Due to a multitude 
of scales, in Sec. 4 we highlight some constraints on the choice of the parameters, which will facilitate the discussions 
of entanglement measure. In Sec. 5 we introduce negativity as a quantifiable entanglement measure for the symmetric 
Gaussian state, and overview its qualitative behaviors in the current configuration. In Sec. 6, we give a quantitative 
analysis of late time entanglement of the coupled oscillators in the common bath. Finally we summarize our findings 
and their physical implications in Sec. 7. Details of key calculations of the covariance matrix elements are shown in 
the Appendix [X} 


II. TWO DETECTORS AT FINITE SEPARATION IN A COMMON FIELD ENVIRONMENT 


We consider a generic setup for the investigation of quantum entanglement between continuous variables directly 
coupled with each other while both interacting with the same environment. Consider two detectors iS'i _2 located 
at X = ±£/2 which remain at rest throughout their dynamics. The internal degrees of freedom (idf) of S' 1.2 are 
represented by harmonic oscillators xi _2 respectively. In a realistic setting the detectors may represent (harmonic) 
atoms whose electronic energy levels are their idfs. We assume the idfs xi ,2 interact directly with each other with a 
bilinear coupling. In addition, each idf interacts with a common environmental quantum field (j) at the locations of 
the detectors. The quantum field provides an indirect interaction between X\,2 which is non-Markovian in nature, in 
the following sense: Any activity in the idf of one detector will generate a perturbance of the field which is picked 
up by the idf of the other detector. Because the propagation of signals in a field cannot be instantaneous, the finite 
time required engenders retardation effects and the evolution histories of the detectors are intertwined imbued with 
memory. This is what gives rise to the non-Markovianity in the evolution of the idf of the system of two detectors, as 
can be seen more clearly in the same setup without direct coupling as from the results of |21j . 
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Therefore in this conhguration three different types of correlations need be considered in the two detectors-quantum 
field system. The first type of correlations comes from the direct interaction between the oscillators’ internal degrees of 
freedom, and is not directly related to the environment. The other two types are directly related to the environment: 
one results from the preexisting correlations of the environment field at locations of the detectors and the other is due 
to the retarded propagation of the field as a consequence of its interaction with the oscillators. 


A. Configuration 


The action that describes the present configuration is given by 

S[x,4']=Jds - Jdsmaxiis)x 2 is)+Jd^x j{x)(j){x)+Jd^x^d^(t>d^^(j), (1) 

with X = (t, x) the spacetime coordinate. The internal degree of freedom Xi of the detector (i = 1, 2) is modeled as 
a harmonic oscillator with equal mass m and equal (bare) oscillating frequency Wb; cr is the coupling strength for the 
direct coupling between Xi and e is the interaction strength between Xi and 4>. Thus the only feature that distinguishes 
between the two detectors is their spatial locations . Assuming a finite separation between them enables us to see 
the full dynamics of environment-induced entanglement which is absent in most prior work on continuous variable 
dynamics (with the exception of |21jl. The current j{x) coupled to the quantum field (j) takes the form 

2 

j{x) = e'^Xi{t)d^^H^-2■^{t)] (2) 

i=l 

which says that the internal degree of freedom Xi of tho detector Si is coupled with the environment held (j) at the 
location z* of Si. 

We express the evolution of the density matrix of the whole system (two detectors-|-the environmental held) using 
the frameworks of closed-time-path integral formalism. After carrying out the path integrals of the held, we arrive at 
a description of the held-inhuenced two-detector system, governed by the effective action 

S'cg[R,A] = fds^l [r'^(s) • A(s) -fA^(s) - 11 ( 3 )] -I- [r'^(s) • • A(s)-f A^(s) • • R(s) 


+ e‘^ I dsds' \ A'^(s) • Gii(s, s') ■ R(s') -I- ^ A'^(s) • Gnis, s') ■ A(s') )■ , 


(3) 


where the row vector = (xi^\x 2 ^^) denotes the internal degrees of freedom of the detectors in the ± branches of 
the closed time path. The superscript T represents the matrix transposition. Moreover 


X± = R± 2^’ 


n = 


2 

(J Ujt 


- fG{2i,s;zi,s') G(zi,s;z2,s') 
^ ’ '' I^G(z2,s;zi,s') G(z2,s;z2,s') 


(4) 


with two kinds of correlation functions G{x;x') of the held dehned by 

Gn{x]x') = i9{t - t')TT^i^pfi[(l){x),(l){x')]'^ , Gh{x;x') = ^ Tr 0 |p,g{(^(a;), (/)(a;')}| ■ (5) 


The former describes the dissipation effects due to the environment in the evolution of the internal degrees of freedom 
in the forms of a local self-force and a nonlocal inhuence, while the latter summarizes the ramihcation of quantum 
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fluctuations of the environment upon the internal degree of freedom. The trace is taken over the degrees of freedom of 
the environment field, and the density matrix in the trace indicates that the initial state of the field is assumed to be 
in a thermal state at temperature ■ When the environmental state possesses translation-invariant symmetry, the 
correlation function G{x;x') has the property G{x\x') = G{x — x') Or, written in explicit space and time variables 
G(x, t; x', t') = G(x — x', t — t'). The Hadamard function Gh is symmetric in its arguments: Gh(x', x') = Gh{x'] x). 

We observe that in the effective action ^CG we have R'^(s)-A(s) = A'^(s)-R(s) and R^(s)-J2-A(s) = A^(s)-J2-R(s). 
Furthermore its imaginary part can be written as a stochastic ensemble average of a Gaussian distribution by the 
well-known Feynman-Vernon [23j path integral identity. 


exp 


— ^ J dsds'A'^{s) ■ Gh{s, s') ■ A{s') = j I?^7^[^]exp i J ds A^(s) ■ ^{s) 


( 6 ) 


where the stochastic Gaussian distribution 7^[^] takes the form 


V[^] = Af exp 


dsds' ^^(s) 


Gj,\s,s')-as') 


(7) 


and Af is the normalization factor. By doing so, we may rewrite the (complex-valued) coarse-grained effective action 
•ScgIRiA] [37] as a real stochastic effective action [3H| 5'e//[R, A,^], weighted by the distribution The new 

variable ^ satisfies the statistics 


( 0 = 0 , 


m^^{s'))=e^ Gh{s,s'). 


( 8 ) 


Since its distribution is Gaussian, all higher moments can be expanded by the second moment via the Wick’s theorem. 
It can be interpreted as the stochastic driving force on the motion of the internal degrees of freedom of the detectors 
in their equations of motion obtained by varying S'e//, 

SSeff 


= 0 , 


mR(s)-I-mfl • R(s) — / ds'Gi^(s, s') R(s') = ^(s). 


( 5 A(s) 

Writing the equation by components may better reveal their physics contents: 

mxi{s)+mu}'^Xi{s) + max2{s)-ds' Gfl(zi, s; Zi, s')xi(s')-k Gi?(zi, s; Z2, s')x2(s') = 0 (s), 

Jo 

mx2{s)+mu:lx2{s) + maxi{s)-e^ I ds Gfl(z2, s; Z2, s')x2(s') + G'i?(z2, s; Zi, s')xi(s') = 6 (s)- 


(9) 


( 10 ) 


( 11 ) 


This set of equations of motion (101 and (11) is of great physical interest in that it contains 1) local external harmonic 
potential, 2) direct interaction between the internal degrees of freedom of the detectors, 3) self-force on each detector 
4) stochastic driving force and 5) indirect, retarded interaction between the internal degrees of freedom via the field. 
The last three effects result from the idf’s coupling with the environment, and in particular, the last one where the 
two detectors idf interact indirectly, as mediated by the environment, which is totally different in nature from their 
direct coupling. This retarded interaction term is weaker as two oscillators are further separated apart. Note that if 
the oscillators carry a charge and are coupled with a quantized electromagnetic field, the corresponding retardation 
interaction term also accounts for the effects of the electromagnetic field produced by the nontrivial motions of the 
other charged oscillators. 
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1. Properties of the Retarded Green’s Function 


We first examine the retarded Green’s functions Gii{x;x') of the scalar quantum field 4>{x), as they describe the 
nonlocal retarded influence of the environment. 


Gr{x; x') =i0{t- t') [4>{x),4>{x')] = ^ 9 {t) 5 {t - r) - Sir + r) 


( 12 ) 


with T = t — t' and r = |x —x'|. We first look into how Gii{x;x') propagates between Xi- Eq- (12) implies the nonlocal 

e{s - i) 


expressions in, say, ( 10 ) will have the effects, 


J ds' Gfl(zi,s;z 2 ,s')x 2 (s') = ^ y ds' 5{s - s' - tj X 2 {s') = X 2 {s - 1) 


(13) 


with = |zi — Z 2 I. The theta function ensures causality, namely, starting the motion at s = 0, then xi will not receive 
the retarded influence from X 2 until at least s = 1. On the other hand in the limit r —i 0, 


Gjt(x,t;x,t') = - —0(t)S'(t), 


and thus (14) in part will account for a purely local effect, only at the locations of the detectors, 

/ ds' Gji(zi,s; zi, s')xi(s') = 

Jo 


1 

27r 


d{s - s')Xi(s') 


- f ds'S{s-s') x[is') 
'=0 Jo 


(14) 


(15) 


due to the reaction of the radiation field. The first term on the righthand side is not well-defined unless we introduce a 
energy cutoff for the environment, whereby the resulting cutoff-dependent expression can be absorbed as a frequency 
normalization of xi- In other words, the natural frequency of x undergoes a shift due to its interaction with the 


environment. Now (15) becomes 

J ds' Gi?(zi,s;zi,s')xi(s') = - ^^2 ’ 

Putting them back into the equations of motion, we arrive at 

Xi(s) -b 27Xi(s) - 27 ^^^^-^X2(s -^) -b a;2xi(s) + crX2(s) = ^ ^i('S). 
X 2 (s) + 27x2(5) - 27 ^xi(s - £)+uj‘^X 2 {s) + crxi(s) = — 6 (s), 

t m 


(16) 

(17) 

(18) 


where 7 = e^/STrm is the damping constant, and uP = tol -\- is the renormalized natural frequency. Formally this 
frequency shift is written as Suj^ = —4xS(0), but we can put it in a form containing the cutoff scale A of the 
model in consideration. 


..^ = -47.(0) = -^ = r 

^ J —00 ^ Jo 


~ 47 47 A 

dK = -/ dK = - 

TT Jn TT 


(19) 


S. Idf Variables y±.- Symmetric (Sum or Center of Mass) and Antisymmetric (Difference or Relative) Normal Modes 

We can decouple this set of equations by superpositions X-i- = (xi + X2)/2 and X- = Xi ~ X2j 

W + 27 X+(s) - 27 x-s(s -£)+ujI X+(s) = — f-h(s), 

t m 


(20) 
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X- +27X-(g)+ 27 X-(s) = 


( 21 ) 


with = (^1 + ^ 2 ) 72 , ^_ = — ^2 and u)\ = ± cr. In the small-separation limit io±l <C 1, we Taylor-expand the 

retarded expressions into x±(s — €) — Xih)^) ~ ^X±('r) + ■ • • and observe that 

X±(s) + 27 X±(s) -h w^X±(s) + 27 ^ x±(s - £) 

27 r ■■■ 

^ X±(s) + 27 X±(s) + (Wfci +( 5 a; 2 )x±(s) + -j[x±{s) - ix±{s) + y X±(s) - X±(s) + • • 


- (1 + 7^) X±(s) + 27(1 ± 1 ) x±(s) + ro±x±(s) ± ^ X±(s) H-, 


4^ A. S'y 

with ,-=F -+ and , = a;^ ± cr. 

7 r 


( 22 ) 


3. Severely Restricted Validity Range of Zero-Separation Results 

In choosing At! = 7r/2, we obtain the equations of motion found in another often used configuration for the investi¬ 
gation of the thermal entanglement in a common bath nanD] 


X+ + ^lX+{s)+ ujI X+is) = — Ms), 

m 

X- + M X-{s) = — C-(s), 


(23) 

(24) 


where ut = ojI_ but = ujI_^ — In that configuration, two coupled oscillators are assumed to be placed at 

the same spatial location and interact with a common thermal bath. With the retardation effects completely ignored, 
the relative mode becomes a undamped driven oscillator. This allows for the initial information of the system to 
persist and its effects are oscillatory in time. This feature fundamentally changes the characteristics of the late time 
entanglement in a common bath, in comparison with the finite separation case we are studying here. 

We can derive the validity range of prior results of the same-location configuration. Keeping the third-order time 


derivative in the small-distance expansion (22), we find the equation of motion for the relative mode becomes 


(l + 7^)x-(s) - ^ X-(s) + (wT. + y)x±(s) + ••• = ^^-(s) ■ 


(25) 


The third-order time derivative term now plays the role of a friction force. Compared with the corresponding equation 


of motion for the CoM mode, the damping of the relative mode is much weaker by an order of (wf)^. Based on (23) 


and (24) we can say that the results for the configuration of two oscillators in the same location interacting with a 


common bath can at best be valid in the very short time before the dynamics of the relative mode is fully relaxed, 
that is, within the time scale much shorter than (ywT.^^)”^. 

B. Solutions and Dynamics 


We now investigate the dynamics of the system derived from the solutions of the equations of motion (10) and (11). 
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1. An iterative scheme to obtain early time solutions 


As a brief interlude, we first mention a method to obtain early time solutions. This gives a clear depiction of how 
mutual influences are transmitted back and forth between the internal degrees of freedom of both detectors. 

The equations of motion of the normal modes are second-order differential equations with fixed delay. One way to 
solve the equations of this type is to make use of iteration. We divide the evolution time t into intervals of length i, 
that is, t £ [{n— 1)1, nl] with n € N. When t falls in the first interval t £ [0, £], the delayed term in (10) and 0 has 
no effect, so the solutions to (101 and 0 can be found exactly. When time evolves to the next interval t £ [£, 2 £] , the 
delayed term has to be taken into consideration, but since t — i£ [0, £], we can substitute the previously found solution 
to the delay term. In so doing the equations of motion in the time t £ [£, 2 £] become ordinary differential equations 
and they are exactly solvable. We may proceed with the same procedures to the next time interval t £ [2£, 3£] and 
further on. This technique is straightforward but becomes very cumbersome after a few iterations. Besides, in this 
iteration scheme, it is not straightforward to identify the existence of the relaxed motion, so it is merely suitable for 
finding the early-time solutions to (10) and 0. Since we are interested in the late time behaviors of the solutions, 
we will not pursue this approach any further. 


2. Systematic method for seeking the solutions 


To systematically find the solutions to the equations of motion (10) and 0, we may perform the Laplace trans¬ 
formation on the equations of motion 


-X± (0) - s x± (0) -k x± (s) 


■27 


-X±(0) -k sx±(s) -ka;^x±(s)T X±(s) = — ^±(s), 


-k 27 s -k a;± T ^ e x±(s) = — C±(s) + [x±( 0 ) -k (s -k 27 ) x±( 0 ) 
where we have used the following identities 

pOO 

/ dt x{t) = -x( 0 ) -k s x(s), 

^0 

poo 

/ dt x{t) = -x( 0 ) - s x( 0 ) -k x(s). 
^0 

poo 

/ dte{t-£)x{t-£)e~’'* = e~"^x{s), 

Jo 

and the Laplace transformation of the function x{^) is defined for t > 0 hj 

pOO 

X(s) = / df xWe"'**-, 

^0 


Solving (261 for x±(s) yields 


X±is) = x±i0) + x±i0) + — ^±is), 


(26) 

(27) 

(28) 

(29) 

(30) 

(31) 


where 




di^\s) = 


9±is) ’ 


(32) 
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and g±{s) = + 27 s + ^ e From the dependence on the initial conditions, we identify that the first two 


terms on the righthand side of (311 correspond to the homogenous solutions to the equation of motion while the third 


term correspond to the inhomogeneous solution. We can find the time evolution of x± (t) by performing the inverse 


Laplace transformation on (31), 


X±{t) = ^ J^ds x±(s) e'**, 


(33) 


where the Bromwich contour F is so chosen to make the integral in (33) well-defined. 

|. All we need to do is to identify the poles 


We use the residue theorem to evaluate the contour integral in 
associated with d) 2 i^)- This is equivalent to finding the solutions to g±{s) = 0. Compared with the corresponding 
equations in the private bath case, the function g±{s) has an extra term proportional to and this renders g±{s) 
a transcendental function, whose exact analytical solution is typically hard to come by. Some approximation schemes 
are needed. 


1. Strong Damping Case When 2ji ^ 1 and 7 > w, the function g±{s) = 0 reduces to 

+ 2x£y + T 27^ = 0 —^ 2xiy + T 27^ = 0 , 

because \y\ must be smaller than 27 ^ otherwise g±(s) will be too positive to be equal to zero. Solving 

27 £?/-ba;±£^T 27 £e"^ = 0 


(34) 


gives 


y = -^ + ^{±e , 

27 


(35) 


where 'I'(z) is the principal solution for x in z = xe^. If we look for a real solution then the relative mode does not 
have any, because 4'(z) is a complex function of z when z < 0. On the other hand since 4'(z) is real for z > 0, it 
is straightforward to see that y becomes positive when 2j > uj^£. A positive real solution to g(s)± = 0 denotes a 
runaway solution when we apply the residue theorem to the inverse Laplace transformation like (|33|). The culprit for 


the existence of the runaway solution can be seen explicitly to be the delayed term in the equations of motion ( 10 ) 


and (11 1 . Thus we can highlight a stark contrast between the Markovian and the non-Markovian motion. 


Markovian dynamics For the Markovian motion of the oscillator, the function g±{s) typically takes the form 

5 ±(s) = -|- 27 s -I- , the Markovian case, (36) 

so the strong damping 7 > uj± merely causes overdamped motion of the oscillator since in this case 5 ±(s) = 0 gives 

(37) 


s = -j±. <0. 
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non-Markovian dynamics However, in the non-Markovian case as shown previously, strong damping can imply a 
real and positive s and in turn unstable motion. In principle, a large damping constant 7 is supposed to efficiently 
damp out the oscillator’s motion. On the other hand, it also means that the oscillator interacts strongly with the 
environment. Thus when two oscillators couple with a common environment, their mutual influences will last and the 
exchange between one another goes on much longer in the course of time evolution. On top of that, when 


uj'^i 


each reciprocal mutual influence is enhanced by a factor —^ > 1 and is accumulated in such a way that it adds 


up constructively, especially for the CoM mode. This acts counter to the damping due to the dissipative force and 
finally wins to become a runaway solution. In contrast, the Markovian motion does not have the retarded terms, so 
the corresponding motion just rapidly damps away in the strong damping case. 


2. Weak Damping Case We first locate the solutions to g±(s) = 0 for large Islf » co±i, and write g±{s) as 




27 - 

1 + — + ^ 
S J 




y-^ev ~ ± 27 ^ = 27£e*(2"+5T|)^ 


with y = si and n S Z. Taking the logarithm of both sides of (38) gives 


2 / + 2 In ~ In 2'yi + * 7 r( 2 n+ - =F -) . 


Now substitute y = u + iv = e® ® with u, n, 0 G K into (39), and we obtain 

u + ln(M^ + v^) — ln 27 £ 


+1 


+ 2 tan-+ X T x) 

u 2 2'. 


= 0 . 


To fulfill (40), the imaginary part of (40) implies Ir^l » |m| such that 


(38) 


(39) 


(40) 


and then it gives v, 


ln(u^ + ~ 21 n|n| 


11 

2n + - T 2 - sgn(u) 


tan ^ ~ sgn(z;) - , 

u 2 


with 


n > 1, 


If we put the solved v back into the real part of (40), we obtain 


u ~ ln 27 £ — 21 n|?;| = In 27 ^ — 21 n( 2 n+^=Fx — sgn(n) 


(41) 


(42) 


(43) 


Unless 7 ^ 1, typically u < 0 for large n. Since the negative values of u will contribute to the decaying behavior of 

(32) once we apply the residue theorem to (33), the factor behaves like 


exp 


u 

-t 


2t . 


- exp<{ -y In- 


2 n + 1 =F 5 - sgii{v) 




TT'I I 


(44) 


Therefore at late times t ^ i, the contributions to the solution (33) from the large n poles are heavily suppressed. 


Alternative, we may compare \u\/i with the damping constant 7 . From the argument that follows (43), we see that 
we typically have |u| yf, in particular for the case of weak damping 7 / 0 ; <C 1. We conclude that the factor will 
fall off much faster than for large n at late time. 


When both oscillators are far away from one another, the delayed term in the equations of motion (10) and (11) 


is suppressed by the factor 7 /.^, so we may treat this delayed term as a small perturbation. Note that in this case 
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since the non-Markovian or the retardation effects are very feeble, the configuration and the equations of motion 
are almost identical to those in the case of the private baths of the same temperature. However, the similarities are 
merely superficial because there is no environment correlation in the private bath case. 


3. Late time behavior At this point we see that the dominant late-time contribution to (33) should come from 
the pole with the smallest value of |u|. We first consider the case 7 <C which covers the limits of weak damping 
and/or large distance between oscillators. The delayed term will be treated as a small perturbation. Let the function 
g±{s) be written as 


~ ^ _ 2/ , 2 , 27 

g±[s) - ^ + + — 


2 /Te 


(45) 


with y = si. The zeros of g±{s) to the lowest order in ^ji is simply given by iu:±i and the minus of them. Here we 
only demonstrate the iui±i case, and the result for the minus sign case can be found accordingly. Next if we take 


into account the contribution of the order 7 /^ in (45), we should expect the correction to the solution is also at least 


of the order 7 /^. Writing y = iui±i + — (a + i/3) + 0( —with a, S K, and substituting this into (45) we can find 
the corrections to the solutions due to the inclusion of the terms of the order 




a = -- 


UJ± L 


uj±i ± sm.uj±i 


P = =F— cosuj±l. 

UJ± 


Thus the fast mode that the zeros of g+{s) are given by 


s+ = ± i w+ + 7 


sintLi+£\ cos u)+i 


UJ+i 


= ± i w+ 


1 - 


7 cosa;+£ 

UJ+ io+£ 


1 + 


sin 

U)+i 


up to the order 0( —)^. Likewise for the slow mode we obtain the zeros of g-{s) with 


s- = ± iw- + 7 


sinw_£\ , cosuJ-t^ 


UJ-k 


= ± i uj- 


1 + 


7 COSW-l 
UJ- LU-l 



^ sinw-^ 

-7 

1 - 7 - 


1 \ 


(46) 


(47) 


(48) 


We observe that the contribution from the retarded term is smaller with large separation i between oscillators, which 
is consistent with our physical intuition. We may group the expressions in sj- in such a way to define the effective 
damping r± and the effective oscillating frequency H±, 


r± = 7 


i± 


sin 


7 cosa;±.^ 

UJ± iij±i 


The effective damping constant is always nonnegative and can only for the relative mode in the limit uj-i 
however such a small separation can induce instability in oscillators’ motion because 

27 COSWi^ 


Wl=uji 


W±=0J± It — • (49) 

-)■ 0 ; 

It—, (50) 

becomes negative when 27 > Negative implies a runaway solution. In addition the condition 27 > uj'^i 

also hints at the possibility that the non-Markovian motion tends to be unstable for strong damping. Therefore in 
the following discussion, we will restrict ourselves to the condition 27 < uj ±£ so that we have a well-defined effective 
oscillating frequency. Finally we remark that this condition is relatively more stringent for the relative mode because 
when the inter-oscillator coupling a is strong, the relative mode frequency uj- tends to be small and it makes the 


condition harder to meet. Furthermore, the perturbative descriptions in (49), though simple in form, tend to generate 


spurious wiggling with £, due to the sinusoidal expression in (50). This wiggling will diminish with large values of £. 
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3. A third method 


Yet another useful approach in dealing with the non-Markovian dynamics is to directly Taylor-expand l/g±(s) 

27 

in terms of the small parameter - and evaluate the corresponding integrals, instead of resorting to the residue 

uj±£ 


. . 27 

theorem. This is possible only when we assume weak coupling with —^ < 1. The idea is to write l/g±{s) as 


1 

5±(s) 


1 


+ 27 s + =F 


27 


■E 


+ 27s + wi ^ + 27s + wi 

^ n—0 ^ 




Owing to the fact 


1 

f ly 1 

■ 1 d ' 

n ^ 

(s^ -I- 27 s -I- 0 ;^)"+^ 

^ 2 ) n\ 

uj± duj± 

s^ -h 27 s - 1 - 


(51) 


(52) 


we find (51) will take the form 


^ 00 

V 


n ^—nsl 

■ 1 d ' 

n ^ 


n\ 

uj± du}± 

s^ -h 2js -\- uj± 


(53) 

The calculations based on this expansion seems free from the artifact wiggling. This approach in a sense is similar 
but more general than the perturbative approach because the perturbative approach can be viewed as a convenient 
case that s in is replace by —iuj±. 

For weak bath-oscillator interaction, the mutual non-Markovian influence between oscillators decays quickly after 
a couple of exchanges, so both approaches will work nicely for the late-time dynamics. This is also the reason why 
at late time the non-Markovian dynamics can be equivalently summarized into a seemingly Markovian form with the 


effective damping and the effective frequency given by (49). However, we would like to stress that the bath in the 


common bath case is correlated at the locations of the oscillators, while the baths in the private bath case are not. 
This is reflected in the separation dependence of d!f ^2 (■s) in common bath case by the perturbative approach. 


III. COVARIANCE MATRIX 


From (31), the general solutions for the motions of the CoM and the relative modes are described by 

X± (t) = (t) x± (0) -f 4"^^ (t) x± (0) -f — [ dt' di^'> {t - t') {t '), 

m Jo 


(54) 


with d^^O) = 1, (i4Ho) = = Oj 4'^^(0) = I; ^nd they are all equal to zero for t < 0. Since we only 

consider the case 2"f/u}'^£ < 1 , we have seen that 4^^('^) exponentially decays with t. At late time after the motions 


of the normal modes are fully relaxed, Eq. (54) asymptotically approaches 

rt 


X±{t) = — [ dt'di^\t-t')^±{f), 
m Jo 


(55) 


so the initial information contained in x±(0) and X±(0) will not survive at late time. Eq. (55) enables us to compute 


the late-time values of the covariance matrix elements. Before proceeding, we first express the force-force correlations 
in terms of the Hadamard functions of the environment, 


o 2 r 


(^+(t')e+(t")) = ^([6(i')+6(i')] [6(i")+6(0]) = y[Gff(0,t'-O + Gi/(zi-z2,t'-t") 
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Ut') - Ut') un - un ) = 2e2 GniO^t' - t") - Gh{^i - ^2,t' - t") 


Ut') + Ut') 6(0-6(0 ) = o, 


where we have made use of the facts that = (^i + 6)/2, C+ = 6 — 6 and 


(6(^06(6)) = e^G_f/(zi,t';zi,t") = e^G//(0,t' — t "), 
(6(^06(^")) = e^G_y(z2,t';z2,6) = e^G//(0,t' — t ”), 
(6(t')6(i")) = (6(t')6(6)) = e2Gff(zi,t';z2,t") = eOff(zi - Z2,t' - t"). 


(56) 

(57) 

(58) 


Thus we are ready to set up the building blocks to construct the elements of the covariance matrix. All we need 
are (x±(t)), ^ ({x+(i)jX-(i)}) and the counterparts for the momentum variables. We begin with (x±(t)) and 


(x+(^))‘ We rewrite X+(0 terms and ^2 and arrive at its late-time value 

(xi( 0 ) = A - n (o(i')o(6)) 


dt'dt” GH{0,t' - t") + Gh(zi - Z2,t' - t") 


2to2 

In the limit t —> oo we arrive at 


(59) 


lim_(xl(i)) = 


t—>-oo 


2m? 


^ |d4^(K)p[Gi/(0, k) + Gh(zi - Z2,k) 
, 27r L 


(60) 


where G//(0 ,k) are the Fourier transforms of d!'^\t' — t"), G//(0,t' — t"), for example. We define the 

Fourier transformation of a function /(r) by 


f{K)= dTf{T)e^^\ 


fir) = 


dn 


.2.71 


/(«)' 


(61) 


(X-(^))- In like manner we find (x?.(t)), 


{X-{t)) = \ [ dt'dt" d^^ ^(t-t')4 H^-6) (C-(t')0(6)) 

7o 


2 e 


2 /•* 


= ^ / dt'dt" di~\t')di~\t") GH{0,t' - t") - Gnizi - 22 , t' - t") 




(62) 


so that at late time we arrive at 


Op 2 roo 1 _ _ 

lim(x^_(t)) = 4 / ^4l4 6«)P Gj^( 0 ,«)-Gff(zi-Z 2 ,K) 


t—>-oo 


^ 27r 


(63) 


• {x+{t)x-{t))- Finally we show that {x+{t)x-{t)) vanishes, 

(x+(t)x-(t)) = 4 f dt'dt" di+\t-t')di-\t-t''){^+{t')C.{t''))=0, (64) 

JO 

due to (^+(<')^_(t")) = 0. The same holds for {x-(t)x+it)), that is, (x-(t)x+(^)) = 0- 
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The displacements xi ,2 of both oscillators are related to x± by 

1 


Xi = X+ + 2 X- 


X2 = X+ 


X- 


(65) 


so with (60l, (63) and (64), we can compute {xi{t)), (xKO) ^od - ({xi(^))X 2 (i)}) in the late time limit. At late 


times they are given by 


lim (xi(t)) = lim {x+it)) + (x+(t)x-) + j (X-(t)) 

t—^oo i—>-oo L 4 . 

lim {xl{t)) = lim \{xl{t)) - (x+(t)x-) + j (X-(t)) 


= lim 

t—yoo 

= lim 


{xl{t)) + ^{x^.it)) 
{xlit)) + l{xUt)) 


lim ^ ({xi(t),X 2 (t)}) = lim (xi(t)X 2 (t)) = lim [(x+(t)) - \ (X-W) 


>-oo 

.2 


( 66 ) 

(67) 

( 68 ) 


Since p± = mx±, we can easily find (p^(t)), {p_{t)) and (p+(t)p_(t)) at late time from the corresponding expressions 
for x±, 


p2 pOO 7 __ _ 

lim ( p ^( t )) = ^ / ^« 2|^(+)(^)|2 Gh { 0 ,^)+ Gh { z ,- Z 2 , k ) 

2 _ _ „27r 


lim {p^_{t)) = 2e^ 

t—¥00 

lim {p+{t)p_{t)} =0, 

t—¥oo 


dn 
, 27r 


|4 T«:)p G//(0 ,k) - G//(zi - Z2,«;) 


(69) 

(70) 

(71) 


and then 


lim {p^{t)) = lim (p+(t)) + - {p_{t)) 

t—^oo t—^oo L 4 

lim (plit)) = lim [(p+(t)) + j {p^_{t)) 
t—>-oo i—>-oo L 4 

\ 4foW’P2(t)}) = ^lim (p+(t)) - ^ (p?.(t)) 


(72) 

(73) 

(74) 


At this point it pays for later discussions to take a closer look at \d^ 2 "'^ Gh{0, k) and Gh{zi — 22 , k). 

We first note that d^^^) is related to d^\s) in (32) by setting s = —i k, thus |d4^('^)P nre even functions of k. 
As for Gh(R, k), we start with the expression of Gh(R, t) 

1 


GH(R,r) = 


f°° d^k 

1 

coth ^ 

g+zk-R,— ZKT _j_ ik-R,+iKT 

1-00 (27r)3 

2 k 

2 

L J 


dn 1 sinwi? , Bk 

-coth — , 

, 27r 47r i? 2 


(75) 


from which we deduced that its Fourier transform is given by 


7^ r-D \ 1^ sin kR /3k 

Gff(R,K) = ^coth — 


(76) 


with K = |k| and R = |R|, and that Gi/(R, k) is also an even function of k. In fact, by the fluctuation-dissipation 
theorem 


^ /.r., X 1 sin kR , Bk , Pk ^ ^ , 

Gi/(R, k) = - - - — coth P— = coth P— Im Gi^(R, k) , 

47r HZ Z 


( 77 ) 
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all uncertainty functions of the canonical variables for the CoM/relative mode can be further simplified, 
for example, the uncertainty of y at late time reduces to 


2 


2/ r' r m2^ / ^ dK coth^ImGfi(K) 

{x {oo))=c^ 7T-\d2{i^)\ Gh{k)=c^ — 


m'^ J_^ 2 Tr ' m'^ ^ G_r(k) + 


2c 

m 


. , Bn 

im / — 02 cotn, 

JO ^ 


where c = 1/2 for the CoM mode and c = 2 for the relative mode, and we have used the fact that 

d2M =--■ 

Here G contains both G(0, k) and G(R, k). Likewise, the late-time uncertainty of p can be written as 

(p^(oo)) = 2mc Im 


f°° dn 2-7 r \ 7 l 3 n 

K 02 (k) coth — 


27r 


Eqs. (78) and (80) apply to both the CoM and the relative modes of the oscillators. 


B. |{te(t),P.(t)}) 


As for the covariance matrix elements like ^ {{Xi{t),Pj{t)}), we have 

\ {{xi{t),Piit)}) = {x+{t)P+it)) + ^ {X-it)p+it)) + ^ {x+{t)p-{t)) + i {x-{t)p-{t)) , 

^ {{X2{t),p2it)}) = {x+{t)p+{t)) - i {x-{t)p+{t)) - i {x+{t)p-{t)) + ^ {x-{t)p_{t )), 

\ ({xi(0>P2(0}) = {x+{t)p+{t)) + \ {x-{t)p+it)) - \ {x+{t)p-(t)) - ^ {x-{t)p_{t )), 

\ ({x2(0,Pi(0}) = {x+{t)p+{t)) - \ {x-{t)p+{t)) + \ {x+{t)p-(t)) - ^ {x-{t)p-(t )). 

At late time t —)■ oo we readily find 

{x+it)p+[t)) = — [ dt'dt" - t')S 2 \t - t”) 

m Jo 

2 ft 

= — / dt'dt” Gff( 0 ,t' - t”) + Gff(zi - Z2,t' - t”) 

2m Jo L 


t —>■ 00 


k|4"''V)P GH(0,/t)-f G/f(zi - z2,/v) 


2m j_^ 27r 

Since |d 2 ~^^(K)P, Gi/(0, k) and Gh{^i — Z 2 ,k) are all even functions of k, we conclude that 

{x+{t)p+{t)) = 0. 

Next we examine {x-{t)p-{t)) and find that 

(X- {t)p- (0) = — [ dt'dt” 4”^ (i - t')d^2~^ {t - t”) {i- {t')^- {t")) 

^ JO 


/o 

2 e2 


= — / dt'dt” d!'~\t')jt\t”) G//(0,t' - t") - Gh(zi - Z2,t' - t") 
m Jo L 

pp2 /-oo 1 _ _ 

t OO -)> * - / — K |4 Gi?(0,K) - G//(Zi - Z2 ,k) 

m 27r L 


Generically, 

(78) 

(79) 

(80) 

(81) 

(82) 

(83) 

(84) 

(85) 

( 86 ) 
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= 0. (87) 

Finally we calculate {x+{t)p-{t)) and (x_(t)p+(t)), 

{x+{t)p-{t)) = - fdt'dt'' di+\t - t')d[-\t - t") 

'm Jq 

= 0 , ( 88 ) 

because (^+(t')^_(t")) = 0. Thus all of the cross-correlations between Xi ^^nd pj vanish at late times, and the only 
remaining nonvanishing elements of the covariance matrix in the limit t —>■ oo are {Xi, 2 {t))j 2 ({itiCO? X^it)}), {p\, 2 Wl 

and - {{px{t),p 2 {t)}). These nonvanishing elements of the covariance matrix can be assembled from the expectation 
values of {x±it)) and (p^(t)). Their computations are shown in the Appendix [A| 


IV. PARAMETER CONSTRAINTS 


To relate to realistic physical conditions it is important to ascertain the range of validity of the most relevant 
parameters in this setup, such as the separation £ between the oscillators, the damping constant 7 , and the inter¬ 
oscillator coupling a. 

We have two oscillators residing in a common bath, which is described by a quantum scalar field initially in its 
thermal state, so they will experience different but correlated field strengths. If we move them closer to one another. 


we find in ( 22 ) that when u!±£ <C 1 , their equations of motion correspond to a configuration in which two oscillators 
seem to be placed at the same spatial location. This is a common configuration used for the investigation of entangle¬ 
ment dynamics |18H2n] . However if we move the oscillators closer to one another, instability in their motion can occur. 


IV. 1. Close Proximity Instablity As discussed in the paragraph below (35), the separation t must be greater than 

27 /w^; otherwise the equations of motion will have a runaway solution. As seen from the criterion, this instability 

7 

can also occur for the strong damping case with a moderate separation between oscillators, i.e. — > ui±l. In fact if 
we write the stability criterion explicitly in terms of the coupling constants, 




27 


± cr 


(89) 


we may easily conclude that the minimum separation to avoid the instability increases with stronger damping and 
stronger inter-oscillator coupling. In other words, strong coupling, either between the oscillator and the environment 
or between the oscillators, tends to induce unstable motion of the oscillator in the current configuration. It is also 
obvious that the slow mode places a stricter constraint on the minimum separation to ensure a stable evolution of 
the oscillators’ motion. This is a distinct feature of the non-Markovian dynamics of the detectors in a shared bath 
configuration, in that the strong damping in the Markovian dynamics in the private-bath configuration only results 
in overdamped motion. 

If one is interested in how thermal entanglement can be sustained for a finite separation between the two detectors 
in a common bath, the <C 1 case is not of special interest because the two detectors essentially merge together. 





22 


The more relevant case would be for uj±£ 1. 


IV.2. Dynamics at late times Based on the previous considerations, the perturbative solution to g±{s) = 0 will 
give the dominant contribution at late time and the approximation will be improved with smaller ratio of 7 /In 
this approximation, all the non-Markovian effects between the oscillators are summarized into the effective damping 
r± and the effective natural frequency W±, 


r± = 7 


^ ^ sin 
lo±£ 


W± = u}± 



7 cos uj±£ 
uj± uj±£ 


of the oscillators, whence the oscillators superficially undergo Markovian dynamics. All the information about the 
history or the finite separation is hidden in these effective parameters. This may not be too surprising because the 
small ratio 7 /( 0 ;^^) implies short memory lapse. At each exchange of the mutual influence between the oscillators, 
the non-Markovian effect will be suppressed by that factor; thus it will not last more than a few exchanges. This also 
implies that the initial information hardly survives at late time. That is why the dynamics of the oscillators in these 
limits looks so Markovian. 


V. ENTANGLEMENT MEASURE 

We will use the negativity as the entanglement measure because it offers unambiguous quantification of entanglement 
for a symmetric two-mode Gaussian state. It is defined in terms of the symplectic eigenvalue ? 7 < of the partially- 
transposed covariance matrix 

J\f{p) = niax{0, , Em{p) = max{0, - ln277<} , (90) 

zr]< 

where ? 7 < is the smaller of the pair of symplectic eigenvalues (? 7 >, ? 7 <) of the partially-transposed covariance matrix 
Entanglement occurs when < 1/2, and the degree of entanglement is described by the negativity. 


A. Negativity in terms of symplectic eigenvalue ?7< 


The covariance matrix V of a two-mode state is define by 



Tr 


p{X,X^} 


(91) 


where X^ = (xi,Pi,X 2 ,P 2 ) and p is the density matrix of the two-mode state. Thus the covariance matrix of a 
two-mode state is a 4 x 4 matrix, which consists of the uncertainties or the cross-correlations between the canonical 
variables of the state. The separability criterion or the entanglement measure based on the covariance matrix is 
particularly convenient and calculable for Gaussian continuous variable systems because its covariance matrix has 
a finite dimension, in comparison with the infinite-dimensional density matrix, which is commonly used in discrete 
systems. The partial transpose of a covariance matrix can be simply constructed by changing the sign of one of 
the momentum variable in the original covariance matrix. The symplectic eigenvalues of the partially-transposed 
covariance matrix can thus be found by taking the absolute value of the ordinary eigenvalues of the matrix 



i J® • , where 
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J = 


0 +1 

-1 0 


If we write the covariance matrix in a block form 


V-( ^ 


(92) 


(93) 


then the symplectic eigenvalues 77 ^ of the partially-transposed covariance matrix can be expressed in terms of 
these 2x2 matrices A, B, C as 


= 


det A -I- det B 


- det C ± 


' ^ det A -I- det B 


-detC -detv 


(94) 


where alternatively det V can be written as det A det B -f (det C)^ — Tr{A • J • C • J • B • J • C'^ ■ J}. 

In the case that A = B and the matrix C is diagonal, the symplectic eigenvalue r]^ of takes a particularly neat 
form 


V< = V 11 V 22 - V 13 V 24 - IV 22 V 13 - V 11 V 24 I = 


_ / (Vii - V 13 ) (V 22 + V 24 ) , V 22 V 13 > V 11 V 24 , 


(Vii + V 13 ) (V 22 — V 24 ) , V 22 V 13 < V 11 V 24 , 


(95) 


with 

Vii = (x?(i)), 

We readily see that 


V 22 = (pm 


^13 = 2 ({xlW-X2(t)}) 


V24 = 2 ■ 


k’li T k’la — 2 ({xiiXi TX 2 }) — 


V22 ± V24 = 2 {{ pi,Pi ±4*2}) = 


2{xl), +, 

+, 


(96) 


(97) 


are associated with the dynamics of the normal modes of the joint system. Thus depending on the sign of V 22 V 13 — 


V 11 V 24 the symplectic eigenvalue ??< will take different forms 



„2 _ 1 

\{x^-)i<pI), 

V 22 V 13 > V 11 V 24 , 

Likewise we find 

P<-\ 

[{xl){p^-), 

k’22k’l3 < k’llV24 • 


2 J 

\{xl){p^-), 

V 22 V 13 > V 11 V 24 , 

In fact we observe that 


[{x^-){pI), 

V 22 V 13 < V 11 V 24 • 


1 r 


V22V13 — V11V24 = 2 


ixDipl) - (x"-)(4) 


so (98 1 can be summarized into 


(98) 


(99) 


( 100 ) 


ryl = min{(x+)(p?.), {X-){p\)] , 


pI = max{(x+)(p?.), (X-)(P+)} • (101) 


This distinction can be handy in interpreting entanglement. Furthermore, in either case, for the entanglement to 
exist, we would like to have the uncertainties of the corresponding canonical variables to be as small as possible such 
that 77 I can be smaller than 1/4. 
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FIG. 1: Behavior of rj^ as a function of separation i between the two oscillators. The symplectic eigenvalue ? 7 < of 
the partially transposed covariance matrix can be used to quantify the entanglement of a symmetric Gaussian 
state, described by V. The value of 77 ^ is given by the minimum between {x+){p^), and {x^){p+)- When 77 ^ is less 
than 1/4, the state is entangled; otherwise it is separable. There may be two different critical separation where the 
entanglement disappears. The point where the curves of {x+){p-), and {x-){p%) intersect always corresponds to a 
separable state. 


B. Interpretation of Entanglement Measure in Terms of Effective Frequencies 

Since the symplectic eigenvalue is constructed out of (x±) and (p±), let us work out their general behaviors 
with respect to the parameters at hand. Recall the equations of motion for the faster/slow modes are given by 

X± +27X±(s) T 27 X±{s - e) +uj‘lx±{s) = — C±(s). (102) 

I m 

From their frequency representations, 

[-K^ - 7 27 K+± cr) =F ^ e*'='^]x±('«) =—?±(k) , <t 4 > x±(«) = — : (103) 

£ m m 

we may attempt to break the into its real and imaginary parts, and group them with the oscillating frequency 
and the damping constant respectively to form the effective oscillating frequency 

1T± =a;^±crT y cosK^, (104) 

and the effective damping constants 

r± = 7 ± ^ sin kI! . (105) 

Kl 

In so doing heuristically speaking we have two effective uncoupled, damped, driven oscillators. For weak oscillator- 
bath coupling 7 , the curve of (k) has a very narrow resonance peak at about k = uj± . Thus we may replace n 
in Fj- and with its typical value, say a;±. However, from hindsight we find this choice is not adequate for W±. 
Gompared with the result of (x±(oo)) at zero temperature, the better option of the effective oscillating frequency will 
be 

£ ’ 


W± = ± cr =F 


(106) 
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in order to avoid wiggling artifact due to the cos expression. Now based on the viewpoint of the effective oscillators, 
we expect that when the temperature is not too high, the uncertainty of the corresponding canonical variables should 
be about the order (apart from the mass scale) 

{pl)-0{W±), (107) 

and on the other hand at high temperature their values are dominated by the temperature, leading to 

(X^)-Cl(^), {pI)-0[-^). (108) 

Let us for the moment forget about any intricacy due to the cutoff scale because it does not really depend on any 
other parameters in the setup. 

Case W+ > W- With these understandings, when and the temperature of the bath is not high, we 

have 


(X^)<(X-), 


and 


(p-}<(pI}, 


(xI}(p^-}<(x^-}(pI}- 


(109) 


From ( 1011 the smaller symplectic eigenvalue ?;< of the covariance matrix in this case will be given by 
?7< = (x+)(P-) ^ C>(^), whenW'+>W_. 


( 110 ) 


With the help of the definition of the effective frequencies (1061, we see that for a fixed w the effective frequencies 


tend to take the opposite trends, that is, when gets bigger, bF_ becomes smaller, and vice versa. It then implies 
that 77^ tends to decrease with larger values of W+ ; thus it is more likely that entanglement can exist for larger W+. 
More precisely speaking, as seen from (106) with a fixed w, the effective oscillating frequency of the CoM mode W+ 


tends to increase with either larger inter-oscillator coupling cr or larger separation i between oscillators. However 
in this case, the effect of the oscillator separation is usually relatively small, overshadowed by the effect due to the 
inter-oscillator coupling, because IF+ > W- already implies that a > The inter-os dilator coupling will play 

the dominant role in sustaining entanglement, if it exists. Therefore stronger inter-oscillator coupling will benefit 
entanglement when both oscillators are far apart. 


Case W- > W+ For the opposite case W- > W+, in the low temperature regime, we have 

(X^)<(X^+), and <(/_), ^ {xl){pl)<{xl){p^-)- (HI) 

The symplectic eigenvalue is then given by 

= (X-) (P+), when W- > W+. (112) 

Following the same arguments, larger values of W_ will allow entanglement to better survive at late time. This 
implies lowering the inter-oscillator coupling strength or shortening the separation between the oscillators can raise 
the value of W-, thus increase the chance that entanglement can exist. In this case W- > W+ leads to cr < 2j/£, so 
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the effect of the inter-oscillator coupling becomes subdominant. The separation between the oscillators plays a more 
important role. The entanglement is more likely to survive for shorter oscillator separation in the case of vanishing 
inter-oscillator coupling. 


Sign Switch in — W- Following the previous arguments and the definition of the effective frequencies, we see 
when a is roughly equal to 2j/£, the difference W+ — W- can switch sign. There the curves of {x+){P-) and {x-){p+) 
will cross over, so the symplectic eigenvalue 77< will change from one form to the other. In addition, the symplectic 
eigenvalues of the covariance matrix can be degenerate in this regime. Since the separability criterion -\-iQ,/2 > 0 
requires 

(113) 


the degeneracy of the symplectic eigenvalues implies 




with the degenerate value r]> = r]^ = rj. 
system always remain in a separable state. 


Eq. (1141 is always greater than or equal to zero, so when a 


(114) 
2x/i the 


High Temperature Regime Finally in the high temperature regime, in general (x±) and {p\) are very large because 


they are all proportional to the temperature, as seen from (108). Hence it is very difficult to maintain either {x+) (P-) 


or {x-){p+) smaller than 1/4. This explains why in the configuration studied here, namely two detectors in a 


common bath, entanglement disappears at high temperature. Eq. (1081 also tells that since the the product {x+){P-) 
or {x^_){p\) is about the order 


pl^{xl){p%)-0{ 




), 


(115) 


the critical temperature of the entanglement will be at most of the order 


P,W± ~ 0{1). 


(116) 


Higher temperature will tend to render the quantity 13W± too small to make greater than 1 /4. Thus entanglement 
may not exist any more at higher temperature. Moreover, it implies that when Wi|_ > VF_ the critical temperature 
is at most about the order W+, while in the case W- > W+, the critical temperature is of the order W-. 

These are the general characteristics of quantum entanglement in a common bath. We now proceed to a more 
detailed quantitative analysis. 


VI. LATE-TIME ENTANGLEMENT ANALYSIS 

The above qualitative description shows the richness of the entanglement behavior at late times (much greater 
than the relaxation time 7 ~^) arising from the interplay between various physical scales. Here we shall examine the 
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zero temperature case and identify more precisely the aforementioned features. The results for the low and high 
temperature cases will be reported in a separate paper [45]. 

Owing to the multitude of length scales in question, we first address the choices of the ranges of length scales involved 
in Sec. |VI A[ We then analyze the critical separation for the regime the direct coupling dominates in Sec. |VI B| and 
for the regime the non-Markovian effects govern in Sec. |VID[ 

A. Identification of Relevant Scales 


We suppose that within the relaxation time of the system, mutual influence between the two oscillators, in terms of 
retardation radiation, has occurred sufficiently many times. It does not take long for weak bath-oscillator couplings. 
This implies 7 “^ > i. Depending on u ^ w^/2, we have 


a > 



(J < 



9 LO 

2 

>a, 


27 47 27 
ujU ’ 

27 47 27 
wU ■ 


(117) 

(118) 


Thus when c > 1, the stability criterion is automatically satisfied if cr < w^/2, but it has to be enforced if cr > a;^/2. 

On the other hand, when c < 1, the stability criterion can be violated when a > 12, but it still needs to be enforced 

if cr < 10 ^ 12. Following our previous discussions, we know that c ~ C>{1) is the region where the “phase transition” of 

occurs and the dominance between the direct and indirect interactions between the oscillators swaps. When the 

inter-oscillator coupling cr is greater than a;^/ 2 , our previous result implies that it is hardly possible to sustain late 

time entanglement for separations shorter than 27 /cr because the instability may spoil everything. 

27 

Even though we require ^ 5 — < 1 for the sake of stability and convergence of perturbative calculation, the choice of 
wit 

uj-i can be ambiguous because there is no restriction on it. For weak oscillator-bath coupling "f/oj- < 1 we still have 
the freedom to choose uj-i ^ 1. We will return to this discussion later. In addition because 7 /w_ < 1, we may not 
bother to distinguish the resonance frequencies of the normal modes from uj± unless confusion exists. 


B. 1 , > 1 : Direct coupling dominates 

We hrst discuss the case c > 1. This is the regime where the direct coupling dominates over the indirect mutual 
influence. The inter-oscillator coupling effectively correlated motions of both oscillators, so entanglement tends to 
survive over finite separation between them. Since there is no restriction about we will discussion the critical 
parameters for two cases: > 1 and < 1 . 


1 . 


When the inter-oscillator coupling is sufficiently strong, that is c > 1, the symplectic eigenvalue ry^, determined by 
(X+( 00 )) (pi ( 00 )), is given by 


2 5(w,7,o-,A) 


1 


-1 7 27 1 87^ 1 /oy ^ '1 


(119) 
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FIG. 2: For sufficiently strong inter-oscillator coupling strength tr, the oscillators may remain entangled even when 
they are far apart. This late time entanglement deteriorates with shorter separation when t > 

Disentanglement will occur when the distance between oscillators is less than the critical separation. We observe a 
sudden transition of at £ ~ In the regime i < if separation is shorter than another critical 

value, the late time entanglement can survive and ameliorate with even smaller separation. This is a consequence of 
the interplay between the direct coupling and the indirect mutual influence between the oscillators. When 
I > 0{2'y ja) direct coupling a dominates, so entanglement in general improves with larger a. On the other hand, 
when (. < 0(2'^jo) the indirect, non-Markovian effect takes over the control of the overall behavior. Finally, we note 
that the curves of r]'^ sink down steeply while moving to smaller values of £ with larger a. 


in the case u!±£ > 1, where < 7 ( 0 ;, 7 , cr, A) is positive and takes the form 

g{u}, 7 , cr. A) = (O?. - 7 ^) cot"^ ^ -f 70 _ In 


A2 




( 120 ) 


with the resonance frequency — 7 ^. The symplectic eigenvalue rj^ approaches a constant when the two 

oscillators are very far apart. If we write this asymptotic constant explicitly in terms of the coupling constants, we 
find 


lim ??< 7 
t—foo 4 


U!^ + a 


+ 7 


'J — a ^ \/w^ + cr 


In- 


A2 


- 1 


+ 0(7"). 


( 121 ) 


27r(a;^ -I- ct) 27r(a;^ -I- cr) 

This constant is smaller with larger cr but grows with increasing 7 because the coefficient of 7 is always positive. It 
implies when the interaction between the oscillators and the bath is not negligibly small, we need sufficiently strong 
inter-oscillator coupling to overcome the vacuum fluctuations of the bath and to maintain entanglement between 
oscillators. For a given strength of the inter-oscillator coupling cr, the damping constant 7 cannot be greater than 


4 1 - 


-|- cr 


7 < 


'JuP- ^ a 

2'k(uP- -|- cr) 


In ■ 




-1 - 


— a 

2'K{bP -I- cr) 


( 122 ) 


otherwise the late time entanglement can not exist. For weaker inter-oscillator coupling ajiip < 1, Eq. (122) is 
approximately given by 

A' 




a; 4(ln A/w — 1) {P 




( 123 ) 
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FIG. 3: The curves of rj^ move upward and to the right with larger 7 , so at a fixed separation this implies smaller 
values of 7 are in favor of entanglement . In the regime ^ > 0{1) the critical separation increases with larger 7 , 
implying that the late time entanglement is harder to sustain. On the other hand, for ^ < 0(1) larger values of 7 in 
fact improves entanglement. Comparing this with Fig. we see the effects of a and 7 on 7 ^, that is, on 
entanglement, are totally opposite in the two regimes c < 0 ( 1 ) and > 0 ( 1 ) . 


This upper bound depends on the cutoff scale A, which originates from the momentum uncertainty. Thus the 
inter-oscillator and the oscillator-bath couplings play competing roles. Intuitively because the damping constant 
7 = e^/ 87 rm is a measure of how much the quantum (vacuum or thermal) fluctuations of the bath can possibly 
disrupt the correlations between the oscillators, so larger 7 tends to make the influence of the bath more destructive 
on the correlations between the oscillators. On the other hand inter-oscillator coupling is expected to increase the 
bond between the two oscillators. Nonetheless, the mechanism of how entanglement is affected by this direct coupling 
is not as intuitive as one would have superficially assumed because when we look into the motion of the normal 
modes of the joint system, we note that it causes the opposite effects. More specifically the stronger inter-oscillator 
coupling will make the CoM mode oscillate faster, which in turn decreases (increases) its position (momentum) 
uncertainty, meanwhile it reduces (amplifies) the position (momentum) uncertainty of the relative mode. However 
from our qualitative discussions based on the effective frequencies, it is not difficult to see how these come into play 
regarding the existence of entanglement. In short, strong inter-oscillator and weaker bath-oscillator interaction benefit 
entanglement only for the case > 1 . 


Since the leading term that depends on the separation £ in (119) is positive, we see that as the separation £ decreases 
the value of rj^ increases as fast as until 7 ^ reaches the value 1/4, where entanglement is about to disappear. 
For weaker but sufficiently strong inter-oscillator coupling, the disentanglement may already happen at comparatively 
large separations. Let £> denote this critical separation, we find 


£> = 


4/3 


l-4a 


(124) 


where a, 13 are functions of w, 7 and cr. 


UJ_ 

a = - -h 

4w+ 


A2 

—LO- + UJ+ (in —2 1) TT^ ( 3 w+ ~ ) + Sw+W- (in —^-l) 




UJ_ 


27rtLi^ 




■7"+ 0(7"): 


(125) 














30 


/3=^7- 


A2 

2 ( 1 - 4 - 1 ) 


0(7^ 


(126) 


with uj^ = ± CT. We will briefly highlight the derivation. In principle the solution of (119) can be found exactly 


because it can be made into a third-order polynomial; however the solutions are too intricate to be of any practical 


use. On observing that the coefficient of the £ term in (119) is relatively small due to higher power in so we 


7 


can drop the £ term for the case uj+£ > 1 and write (119) into the form 


/3 


1 


’'< = “+P = 3 


( 127 ) 


after we Taylor-expand (119) in orders of the small parameter 7 /a;±. Eq. (127) is equivalent to a second-order 


polynomial in £, and one of whose solutions is then (124). At this step, we see the parameters a, P in (124) are still 


quite complicated and may attempt to further expand (124) in terms of small ratio It turns out that it may 

not always work, especially for the case when £ is sufficiently large, the saturated value of the curve happens to 
lie right below the horizontal like 1/4 (please refer to Fig. |^. To meet this situation, the coupling constants 7 and 
a must be so chosen to make the value of for far-apart oscillators close to but smaller than the critical value 1/4, 
that is. 


lim 

n+i>i ^ 4 


(128) 


Since both P and (1 — 4a) are already small, any minute variation may induce large error in approximation. On top 
of that, the 0(1) and 0 ( 7 ) terms in (1 — 4q) are comparable in magnitude but take on opposite signs. This makes 
(1 — 4a) at least roughly of the order 0 ( 7 ^). It may appear odd how the 0(1) and 0 ( 7 ) terms in (1 — 4a) can be 


comparable even for small 7 . But recall that in order to (almost) satisfy (128), the coupling constants 7 and cr must 
take on particular values. The consequence of that is the (almost) cancellation of the 0(1) and 0 ( 7 ) terms in (1 —4a). 
This is not straightforwardly revealed via the Taylor expansion of small 7 and will yield erroneous results for the case 


described by (128). Thus (124) can give a satisfactory estimation of even is large, but it still does not work 


for the case £1±£^ ^ 1 , where the contributions of the order higher than 0 ( 7 ^) must be taken into consideration. 


g. u-e < 1 


As seen in Fig. when the inter-oscillator coupling is stronger, disentanglement could happen at relatively small 
separation. Thus we use the short-distance approximation {ui±£ < 1) of 77 ^ 


ri% = 


UJ_ 

4uj^ 


(u!+ -f a;)a;^ 7 


TTWl 


UJ+ 




UJ2 


+ J 




{uj^ + lo'^)lo_ (w_|_-I-w_)(w^-|-) — w+w/. ( 7 e -I-In A£) 7 


0-0 ^1-1" ' 

a;+ 


40;+ 


7 

J 


7ro;+ 0 ;+ 

_ (cu+-fo;-)-o;+( 7 e + lnA£) 7 
4o;+ 7 ro;+ a;+ ' 0 ;+ 


bo + bi^ 

UJ+ 


7 


7 

Co + Cl-h 


7^ 

^( 4 ) 


7 

oU 


o( 


7 W 1 


+ ^(7772") ’ 


o , 

7^1 


(129) 
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to find the critical separation ^>. In (1291 the shorthand notations are 
,3 0 ,2 I , ,2 


Go = 


w_ 

44 ’ 


ai = 


(a;+ + CLi)tLi^ 


bo — 


Co = 






=- 


(a;+ + oj -) (wl + u!_) — uj+u}_ ( 7 ^ + In A£) 


4tLl4 


Cl = 


(a;+ + w-) — a;+ ( 7 ^ + In A£) 


iu;_|_ TTW^ 

By means of the criterion 7 ^ = 1/4, the critical separation is determined approximately by 




ai 


+ bi 






Cl 


e> = + 


,2 »( 0 ) ' 


2 a( 


a;+ 


0 ~r ^0 


(130) 


where is the leading contribution in the critical separation 




27 


-(w^ + w?.) + [(w2 - )2 + 4 w 3 a;_]l/2 ’ 


(131) 


To arrive at (130) we have made the assumption 7 (jJ± to simplify the result and recall that = oj^ ± a. 


3. Comments 


In this regime the choice of a is limited by two constraints: <; > 1 and 27 < uj^£. They imply that 

27 ^ ^2 27 

-y 


(132) 


For small values of a, the critical length tends to be large and (124) gives a better approximation, while large 


a implies a small critical separation so (130) is more suitable. Within the range (132), the critical separation t. 


decreases with growing inter-oscillator coupling strength. The physical interpretation regarding the separation of the 
oscillators is much less straightforward due to the fact that non-zero separation picks up the correlations existent in 
the bath between two spatial locations and is the cause of non-Markovianity in the propagation of mutual influences. 
Since the non-Markovian influence of one oscillator will propagate by means of the intermediate bath in the form of 
the retarded radiation to the other oscillator, thus inevitable introducing a disparate phase and further degrading 
coherence between them. However, this non-Markovian effect is greatly reduced by the separation. That may offer 
a physical understanding why in this regime entanglement can possibly improve with increasing separation between 
oscillators. 


C. <;~0(1) 

When we further decrease the separation between oscillators, the joint system eventually evolves into a sepa¬ 
rable state. Direct coupling is now not strong enough to sustain the late-time entanglement. We will gradually 
come to another critical separation £x ~ 27 / 0 -, where the crossover between the values of {)^{oo)){p^{oo)) and 
(x?_(oo))(p4oo)) happens. Once the separation is smaller than this critical value, the symplectic eigenvalue will be 
represented by (y?.( 00 ))( 4 ( 00 )) instead of {x\{oo)){p^(oo)) because now the former is the smaller among the two. 
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D. ? < 1: non-Markovian interaction dominates 

In this regime the non-Markovian effects is more important and the effects of direct coupling becomes subdominant, 
so we expect that the results should be similar to the earlier works [mui, where direct coupling is absent. The 
stability criterion implies that 

< £ < — , and a <%- . (133) 

cr 2 

Thus when the inter-oscillator coupling is very weak, i.e., aJuj <C 1, we find can be greater than one. On the other 
hand, if cr/w^ < 1/2, then the value of is very small, at most about the order ©(y/w^). 


1 . 


First considering the case cr/o;^ <C 1, we use the large separation approximation, that is, ^ 1, to express rj^, 


2 _ \ , 
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(134) 


Again it shows that for <^ < 1 but a;_£ ^ 1, the oscillators remain separable until the distance between them is 
shortened to .^ = ©(wZ^), where starts decreasing rapidly. Apparently it is not suitable for identifying the critical 
separation. Hence we need an approximation of 77 ^ that works better in the regime W-l < 1. 

2. uj-i < 1 

We find the small-separation approximation of is roughly given by 


??< - 


oj- u)+ + oj- 7 
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The leading term of the order is negative. In addition, in this regime a < 27 /.^, the stability criterion implies 


7 

ujU 


(135) 


27 27 

-— < - — < 1 

UJ^^ ujU ’ 


(136) 


if we recall uj± = Therefore as far as the leading terms in each pair of brackets in (135) are concerned, we 


find 


1 7 ^ UJ^ +UJ_ -y 


Aoj+oj^ Aoj+oj^ I 


-k w/. 7 
dw+w/. I 


2uj‘^e 


- 1 


-k • • • < 0 , 


(137) 


so 77 ^ decreases with smaller separation i until again it crosses the value 1/4. This means that entanglement is likely 
to exist and furthermore is improved for shorter separation. 
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FIG. 4: By numerically solving for 77 ^ = 1/4, we generate the relation between the critical separation with 
respect to 7 and a for the case described by ( |138 ) in the regime < 1. Late time entanglement indeed improves 
with larger 7 and smaller cr. In this regime, the indirect non-Markovian mutual influence dominates over the effects 
of direct coupling cr between the oscillators. Such history-dependent effects make the physical interpretation less 
straightforward. In the plot we choose the oscillating frequency w = 5 and the cutoff scale 10000. 


Compared with (129) where the cutoff-dependent term takes the form lnA£, the counterpart in (135) has quite a 


‘JJ \ I I 

distinct form. In i. Thus the logarithmic cutoff dependent terms in (135) weigh much more than those in (129) 


unless 7 is vanishingly tiny, and this makes the application of the iteration scheme in (1351 much trickier. On account 


of this difficulty, we observe in (1351 that the contribution in the first pair of brackets on the righthand side can be 


negligible compared with the rest. If we ignore it, the remaining expressions of (135) can be used to formally solve 
ryl = 1 /4, leading to 


= 


(a;3 +a;3)+cc2a;(7, + ln5^^) 


4a;+u; 


TTW+w/. 
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(138) 


Repeated iterations in on the righthand side of (138 1 rapidly produce better improvement on the critical separation 


£<, although the result can become very complicated in a couple of iterations. 


This approximation works best for smaller values of a because the ratio 


27 


will remain small. From the involved 


UJ _ 

expression of ( |138[ ) we can still identify the facts that the critical separation grows with increasing 7 and decreasing a. 
In Fig.j^the results, generated numerically, support this trend. On the contrary, stronger inter-os dilator coupling will 
render this critical separation even shorter, but in the end we will come to a situation that either the approximation 
breaks down or the instability occurs. 

In the limit cr/w^ <C 7 /w, that is, with negligible direct coupling between the oscillators, the critical separation 
approaches a value with very mild dependence on 7 . If the ratio 7 / 0 ; is much smaller than one, then this value 
essentially becomes a constant. 


2uj (In 
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r ^ 7 
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(139) 


as shown in Fig. Since this is the iterative expression, if we truncate it to the first order, we can substitute in 


the logarithm on the righthand side by —. 

2a; 
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FIG. 5: When the inter-oscillator coupling is vanishingly small, the critical separations are almost independent of 
the damping constant 7 . Note that since the inter-oscillator coupling is too weak, there is only one critical 
separation for each curve. It is impossible to sustain entanglement at large separation in this case. 


So far, we note that when we look for the critical separation, we repeatedly come across solving a transcendental 
function of the form 


C2 


-f Cl In £ -I- Co = 0 . 


(140) 


Some of its solutions, as shown previously, can be found by iterative substitutions, but in fact it can be expressed in 
terms of a special function, the product logarithm or Lambert W-function yV{z), which is the principal solution for 


w in z = w e'". Hence for example, Eq. (1391 can be written as 


2uj - 




(141) 


where W(fc, z) is the kth solution to z = we'^ other than the principal one. In addition, the trails of the cutoff scale 


are dotted everywhere and they cannot always been argued away. As in the case of (1391, the effect of the cutoff scale 


persists and is independent of the damping constant 7 even though 7 is already tiny. 

Thus we see, for sufficiently weak inter-oscillator coupling, when the separation between the oscillators is less than 
this critical value, the entanglement between them can survive at late time. If we further reduce the separation, we 
will at last come to the situation that the dynamics of the joint system becomes unstable. 


VII. SUMMARY OF RESULTS 


Here we summarize our analytical results and numerical supports as shown in the figures. We first note that 
different normal modes have different effective relaxation time scales, and these time scales also depend on the 
separation between oscillators. When the mutual separation is small u}±i 1, the motion of the slow mode relaxes 
much slower than that of the fast mode. Thus the late-time results will not be valid unless the motion of all normal 
modes are relaxed. Within the relaxation time, we assume the mutual influence between the oscillators has been 
exchanged for sufficiently many times. Then we find that except for some extreme cases, typically speaking, in the 
regime > 1, the late time entanglement improves with I) stronger inter-oscillator coupling, 2) weaker bath-oscillator 
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coupling, or 3) longer oscillator separation. The critical separation at which disentanglement occurs thus reduces with 
larger value of tr but smaller value of 7 , as can be seen in Figs. [^and[^ In addition, within this regime, the direct 
coupling between oscillators dominates over the indirect, non-Markovian interaction. The latter effect depends on 
oscillators’ histories and is manifested as the retarded radiation mediated by the environment. Thus the magnitude 
and the phase of this non-Markovian effect will depend on the separation. Roughly speaking at large separation, we 


see from the Fourier or Laplace transforms of the normal-mode equations of motion (103) that the retarded term 
tends to introduce the phase difference from the local terms for the contribution of each environment mode due to 
the delay by This out-of-phase effect is factored out by larger oscillator separation, so it more or less explains the 
improvement of entanglement with increase separation in this regime. 

In the neighborhood c ~ C>(l)j comparison of the behavior of between various values of 7 or cr is difficult to 
make because that is the region where the transition of 77 ^ happens. For example, in Fig. when i ~ 0.03, the curve 
of 77 I for cr/w^ = 0.2 is described by (x?.(oo))(p^(oo)), but the curve that corresponds to tr/o;^ = 0.4 is still given by 
(x+(oo))(p^{oo)). Thus we don’t see the comparison between them particularly meaningful. 

When c < 1, we find the entanglement tends to improve with 1) weaker inter-oscillator coupling, 2) stronger 
bath-oscillator interaction, or 3) shorter separation. In other words, stronger bath-oscillator interaction, or weak 
inter-oscillator coupling is preferable for the late time entanglement so that it can still be sustained at larger mutual 
separation. This conclusion can be confirmed in Fig. where the critical separation is numerically found by solving 
=1/4 for different values of 7 and a. These results seem less intuitive because in this regime the non-Markovian 
effect is much more prominent than the direct coupling effect. Except for the case that tr is vanishingly small, the 
typical values of I satisfy a;±£ < 1. This implies the delay introduces only a nominal phase difference between the 
retarded term and local terms for the contribution from each environment mode. This subtle difference can be very 
sensitive to the values of cr and 7 , and then is readily amplified by the tiny value of the mutual separation, as can be 


seen in (103). This points out the elusive aspect of the non-Markovian phenomena. 

However, the effective oscillating frequency IF±, introduced earlier, can be a very attractive tool in deciphering the 


non-Markovian effects. Recall that the corresponding effective frequency W± is defined by (106 1 

IF,J = ± cr T ^ , 


where the non-Markovian effect is already encapsulated in the expression 27 /t'. We see that this non-Markovian effect 
tends to lower (raise) the oscillating frequency of the CoM (relative) mode. The reduction (enhancement) increases 
with larger bath-oscillator coupling 7 and shorter separation i. This can be understood from the equations of motion 
of the normal modes 

X± + 27 X±(s) +u;± x±(s) = ±2j x±(s - £) + — f±(s) • 

£ m 

We see on the righthand side the delayed term in fact is the retarded Coulomb-like force. For the CoM mode, it exerts 
a repulsive force, which counteracts with the restoring force from the harmonic potential of the CoM mode, effectively 
leading to slowdown of the CoM-mode oscillation. Likewise, this retarded Coulomb force is attractive for the relative 
mode, and thus works hand in hand with its restoring force, causing a more rapid oscillation. The effect of the direct 
coupling can also be easily understood by means of this effective interpretation in the same fashion. Earlier on we 
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have argued that when ^ > 0(1), we have W+ > W- and the value of is roughly given by 

27 , 


W-. 




uj — a + 


UJ^ +<T- 


27 


), 


(142) 


if we explicitly spell out the dependence of W± on 7 , a, oj and £. This clearly explains that the entanglement 
will improve with larger inter-oscillator coupling, weak bath-os dilator interaction, and longer separation between 
oscillators. On the other hand, when g < 0(1), we instead have W- > W+. In this case the value of takes 
a different form 


V%-Oi^) = Oi 


2 27 

■-UJ^ + CT-/-| 3 


2 ^ 27 

— CT H—- 


)■ 


(143) 


It also nicely explicates why the interaction and the mutual separation play totally opposite roles in this case, in 
comparison with their effects for the situation > 0 ( 1 ). 

In this paper we have shown that spatial dependence and direct coupling of two oscillators in a common zero- 
temperature quantum field environment engenders very rich features in their entanglement dynamics arising from 
the interplay between the different scales involved. In a sequel paper |45j we will consider the entanglement between 
two coupled oscillators in a common finite temperature bath, and address the question whether and under what 
conditions entanglement can be sustained at high temperatures, as some authors alluded to recently |49j . 
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Appendix A: evaluation of {x±{t)) and {p±(t)) at zero temperature 


Here we first evaluate {x±{t)). In general, at late time the expectation values of lim (y . (t)) and lim {p±{t)) have 

t—^OO t—^OO 


be shown respectively to be 


/ 2 / 2 ci f dK -^(±). , (dK 

(x±(oo)) = -Im / — ^2 (k) coth — , 

m Jq ZTT 2 

(p±(oo)) = 2mc± Im [ — d^^\K) coth — , 

Jo 27r 2 


(Al) 

(A2) 


for the CoM and the relative modes. However, due to the presence of the hyper cotangent function and non- 
Markovianity inherent in d 2 ^^(«;), there is no closed form of the integral for arbitrary values of the parameters 
at hand. Hence we will evaluate (y^(oo)) and (p^(oo)) in the zero, the low- and the high-temperature limit, where 
the factor coth ^ is approximately given by 


1 + 2 e-^«, 



» 1 , 

Pk 1 . 


(A3) 
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In addition, we will assume the coupling between the oscillator and the environment satisfies the condition 27 < up'I 
such that we can treat the non-Markovian term as small perturbation. 

I3 tx 

For the current article, we only discuss the zero temperature case, where coth — takes the value of unity, and we 
Taylor-expand c? 2 ^\k) in terms of the small parameter —^ 


02 («) = 


1 d 


—P — i u!± duj^ 


±4e* 


— P — i 27K + U!^ 


(A4) 


r(±) 


where is obtained if we substitute s = —in into d 2 ^^(s) in ([3^, that is, {k) = {—i k). We observe 

27 

that since in general u!± are different, we will restrict ourselves to the tighter constraint among the condition „ . 
The expectation value (x^(oo)) in the zero temperature limit is given by 
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with c+ = 1/2 for the CoM mode and c_ = 2 for the relative mode. We begin with the calculation of the integral 
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dK 
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where wj — 7 ^ is the resonance frequency for the CoM/relative mode. Next the integral that accounts 

for the non-Markovian contribution can be expressed as differentiation of a simpler integral 
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so we evaluate that integral first 
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Thus we have given by 
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We see that both and 1 ) 2 ^^ are pure imaginary. The exponential integral function Ei(z) is defined according to 


Eiiz) = -I 


00 -S 

ds - 
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The late-time value of {x±(t)) in the zero temperature limit is then given by 


(x±(oo)) = -^Im 
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Similarly the late-time value of the expectation value (p^(t)) in the zero temperature limit takes the form 


, 2 . T f°°dK 

{p±{oo)) = 2mc±Im J — 

We start with the integral 
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which apparently will diverge. Thus we replace the upper limit of the integral by a cutoff frequency A to regularize 
the integral. We find, after regularization, is given by 
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This type of cutoff-dependent expression is often seen in a system interacting with a quantum-field environment, as 
the consequence of an infinite number of degrees of freedom associated with the field. The introduction of the cutoff 
scale is based on the understanding that any effective physical system or model has its range of validity. The cutoff 


frequency thus represents the highest energy scale permissible with the model. The leading term in (A14|, which is 
linear in A, is of no relevance to us because it belongs to the real part of . What is of our concern will be the 
term that is proportional to In A, which has a mild dependence on the cutoff scale, so it does not pose a serious issue. 
On the other hand defined by 
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is independent of the cutoff scale. The integral in .^2 is evaluated to be 
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Again we see is imaginary. Therefore we have (p^(oo)) in the zero temperature limit given by 


(p±(oo)) = ^^Im 
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